import copy
import math

import cv2
import numpy as np
import time


def preprocess(image_path):
    """
    1.彩色加载
    2.转灰度
    :param image_path:
    :return:
    """
    # 读取图片
    image = cv2.imread(image_path)

    # 转成灰度
    image = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)

    return image


def brenner(image_path):
    """
    brenner梯度函数
    :param image:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape
    score = 0
    for i in range(x - 2):
        sub = image[i+2,:]-image[i,:]
        score += cv2.sumElems(sub**2)[0]
    return score / (1.0 * x * y)


def tenengrad(image_path):
    """
    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape

    grad_x = cv2.Sobel(image, -1, 1, 0, ksize=3)
    grad_y = cv2.Sobel(image, -1, 0, 1, ksize=3)

    score = cv2.sumElems(cv2.sqrt((grad_x ** 2 + grad_y ** 2)))[0]
    return score / (1.0 * x * y)


def laplacian(image_path):
    """

    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    x, y = image.shape
    imagelap = cv2.Laplacian(image, cv2.CV_64F)
    score = cv2.sumElems((cv2.absdiff(imagelap,0)))[0]
    return score / (1.0 * x * y)

def tenengrad2(image_path):
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape

    hori = cv2.Sobel(image, cv2.CV_64F, 1, 0)
    vert = cv2.Sobel(image, cv2.CV_64F, 0, 1)

    assement_weight = 1.0
    fx, fy = hori * hori, vert * vert
    dst = np.sqrt(fx + fy)

    dst = cv2.normalize(dst, None, 0, 255, cv2.NORM_MINMAX)

    mean_val = np.mean(dst)

    final_score_ori = np.sqrt(np.mean((dst - mean_val) * (dst - mean_val)))
    final_score = final_score_ori * assement_weight

    return final_score

def SMD(image_path):
    """

    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x,y = image.shape

    kernel_x = np.zeros((3,3),np.float32)
    kernel_x[1,2] = -1.0
    kernel_x[1,1] = 1.0
    kernel_y = np.zeros((3, 3), np.float32)
    kernel_y[0, 1] = -1.0
    kernel_y[1, 1] = 1.0

    smd_x = cv2.filter2D(image,-1,kernel=kernel_x)
    smd_y = cv2.filter2D(image,-1, kernel=kernel_y)

    score = cv2.sumElems(cv2.absdiff(smd_x,0)+cv2.absdiff(smd_y,0))[0]

    return score / (1.0 * x * y)


def SMD2(image_path):
    """
    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape

    kernel_x = np.zeros((3, 3), np.float32)
    kernel_x[1, 2] = -1.0
    kernel_x[1, 1] = 1.0
    kernel_y = np.zeros((3, 3), np.float32)
    kernel_y[1, 1] = 1.0
    kernel_y[2, 1] = -1.0

    smd_x = cv2.filter2D(image, -1, kernel=kernel_x)
    smd_y = cv2.filter2D(image, -1, kernel=kernel_y)

    score = cv2.sumElems(cv2.absdiff(smd_x, 0) * cv2.absdiff(smd_y, 0))[0]

    return score / (1.0 * x * y)

def variance(image_path):
    """

    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape
    mu = cv2.mean(image)[0]
    score = cv2.sumElems((image - mu) ** 2)[0]
    return score / (1.0 * x * y)

def energy(image_path):
    """

    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape

    kernel_x = np.zeros((3, 3), np.float32)
    kernel_x[1, 2] = -1.0
    kernel_x[1, 1] = 1.0
    kernel_y = np.zeros((3, 3), np.float32)
    kernel_y[1, 1] = 1.0
    kernel_y[2, 1] = -1.0

    smd_x = cv2.filter2D(image, -1, kernel=kernel_x)
    smd_y = cv2.filter2D(image, -1, kernel=kernel_y)

    score = cv2.sumElems(smd_x*smd_x + smd_y*smd_y)[0]

    return score / (1.0 * x * y)


def vollath(image_path):
    """

    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    x, y = image.shape
    score = 0

    for i in range(x - 1):
        mul = image[i+1,:] * image[i,:]
        score += cv2.sumElems(mul)[0]

    score = score - x * y * (cv2.mean(image)[0] ** 2)
    return score / (1.0 * x * y)


def entropy(image_path):
    """
    :param image_path:
    :return:
    """
    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path

    numEntries = len(image)  # 样本数
    labelCounts = {}  # 该数据集每个类别的频数
    for featVec in image:  # 对每一行样本
        currentLabel = featVec[-1]  # 该样本的标签
        if currentLabel not in labelCounts.keys(): labelCounts[currentLabel] = 0
        labelCounts[currentLabel] += 1

    score = 0.0
    for key in labelCounts:
        prob = float(labelCounts[key]) / numEntries  # 计算p(xi)
        score -= prob * math.log(prob, 2)  # log base 2
    return score


def EAV(image_path):
    """

    :param image_path:
    :return:
    """

    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)

    x, y = image.shape

    score = 0
    for i in range(1, x - 1):
        for j in range(1, y - 1):
            score = abs(image[i - 1, j - 1] - image[i, j]) * 0.7 + \
                    abs(image[i - 1, j] - image[i, j]) + \
                    abs(image[i - 1, j + 1] - image[i, j]) * 0.7 + \
                    abs(image[i + 1, j - 1] - image[i, j]) * 0.7 + \
                    abs(image[i + 1, j] - image[i, j]) + \
                    abs(image[i + 1, j + 1] - image[i, j]) * 0.7 + \
                    abs(image[i, j - 1] - image[i, j]) + \
                    abs(image[i, j + 1] - image[i, j])

    return score / (1.0 * x * y)

def nrss(image_path):
    """

    :param image_path:
    :return:
    """

    def ssim(i1, i2):
        c1 = 6.5025
        c2 = 58.5225

        I1_2 = i1 * i1
        I2_2 = i2 * i2
        I1_I2 = i1 * i2

        mu1 = cv2.GaussianBlur(i1, (11, 11), 1.5)
        mu2 = cv2.GaussianBlur(i2, (11, 11), 1.5)

        mu1_2 = mu1 * mu1
        mu2_2 = mu2 * mu2
        mu1_mu2 = mu1 * mu2

        sigma1_2 = cv2.GaussianBlur(I1_2, (11, 11), 1.5)
        sigma1_2 -= mu1_2
        sigma2_2 = cv2.GaussianBlur(I2_2, (11, 11), 1.5)
        sigma2_2 -= mu2_2
        sigma12 = cv2.GaussianBlur(I1_I2, (11, 11), 1.5)
        sigma12 -= mu1_mu2

        t1 = 2 * mu1_mu2 + c1
        t2 = 2 * sigma12 + c2

        t3 = t1 * t2

        t1 = mu1_2 + mu2_2 + c1
        t2 = sigma1_2 + sigma2_2 + c2
        t1 = t1 * t2

        ssim_map = cv2.divide(t3, t1)
        mssim = cv2.mean(ssim_map)

        ssim = (mssim[0] + mssim[1] + mssim[2]) / 3

        return ssim

    import warnings
    warnings.filterwarnings("ignore")

    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path

    Ir = cv2.GaussianBlur(image, (7, 7), 6, 6)

    G = cv2.Sobel(image, cv2.CV_32FC1, 1, 1)
    Gr = cv2.Sobel(Ir, cv2.CV_32FC1, 1, 1)

    cols = len(G[0])
    rows = len(G)

    block_cols = cols * 2 / 9
    block_rows = rows * 2 / 9

    best_G = None
    best_Gr = None

    max_stddev = .0
    pos = 0

    for i in range(64):
        leftx = (i % 8) * (block_cols / 2)
        lefty = (i / 8) * (block_rows / 2)
        rightx = leftx + block_cols
        righty = lefty + block_rows

        if leftx < 0:
            leftx = 0
        if lefty < 0:
            lefty = 0
        if rightx > cols:
            rightx = cols - 1
        if righty > rows:
            righty = rows - 1

        leftx = int(leftx)
        lefty = int(lefty)
        rightx = int(rightx)
        righty = int(righty)

        temp = copy.deepcopy(G[leftx:rightx, lefty:righty])

        mean, stddev = cv2.mean(temp), np.std(temp)

        if stddev > max_stddev:
            max_stddev = float(stddev)
            pos = i
            best_G = temp
            best_Gr = copy.deepcopy(Gr[leftx:rightx, lefty:righty])

    score = 1 - ssim(best_G, best_Gr)
    return score

def NR_IQA(image_path):
    """

    :param image_path:
    :return:
    """

    def compute_blur_IQA(gray_image):
        gray_image = copy.deepcopy(gray_image)
        gray_image = np.float32(gray_image)

        kernel_h = np.zeros((3, 3), np.float32)
        kernel_h[0, 1] = -1
        kernel_h[2, 1] = 1
        grad_h = cv2.filter2D(gray_image, cv2.CV_32FC1, kernel_h)

        kernel_v = np.zeros((3, 3), np.float32)
        kernel_v[1, 0] = -1
        kernel_v[1, 2] = 1
        grad_v = cv2.filter2D(gray_image, cv2.CV_32FC1, kernel_v)

        mean = cv2.mean(grad_v)[0]
        mask = grad_h > mean
        mask = mask / 255
        mask = np.float32(mask)
        C_h = cv2.multiply(grad_h, mask)

        edge = np.zeros((len(C_h), len(C_h[0])))
        for i in range(1, len(C_h) - 1):
            for j in range(0, len(C_h[0])):
                if C_h[i - 1, j] < C_h[i, j] and C_h[i + 1, j] < C_h[i, j]:
                    edge[i, j] = 1

        A_h = grad_h / 2
        BR_h = cv2.absdiff(gray_image, A_h)
        BR_h = cv2.divide(BR_h, A_h)

        A_v = grad_v / 2
        BR_v = cv2.absdiff(gray_image, A_v)
        BR_v = cv2.divide(BR_v, A_v)

        int_blur = np.zeros((len(BR_v), len(BR_v[0])), np.float32)
        for i in range(len(int_blur)):
            for j in range(len(int_blur[0])):
                int_blur[i, j] = BR_v[i, j] if BR_v[i, j] > BR_h[i, j] else BR_h[i, j]

        blur = int_blur < 0.1 / 255
        blur = np.float32(blur)
        sum_inv_blur = cv2.countNonZero(int_blur)
        sum_blur = cv2.countNonZero(blur)
        sum_edge = cv2.countNonZero(edge)

        blur_mean = sum_inv_blur / sum_blur
        blur_ratio = sum_blur / sum_edge

        return blur_mean, blur_ratio

    def compute_noise_IQA(gray_img):
        gray_img = copy.deepcopy(gray_img)
        gray_img = np.float32(gray_img)

        gray_img = cv2.blur(gray_img, ksize=(3, 3))

        kernel_h = np.zeros((3, 3))
        kernel_h[0, 1] = -1
        kernel_h[2, 1] = 1
        grad_h = cv2.filter2D(gray_img, cv2.CV_32FC1, kernel_h)

        kernel_v = np.zeros((3, 3))
        kernel_v[1, 0] = -1
        kernel_v[1, 2] = 1
        grad_v = cv2.filter2D(gray_img, cv2.CV_32FC1, kernel_v)

        D_h_mean = cv2.mean(grad_h)[0]
        D_v_mean = cv2.mean(grad_v)[0]

        N_cand = np.zeros((len(gray_img), len(gray_img[0])), dtype=np.float32)
        for i in range(len(gray_img)):
            for j in range(len(gray_img[0])):
                if grad_v[i, j] < D_v_mean and grad_h[i, j] < D_h_mean:
                    N_cand[i, j] = grad_v[i, j] if grad_v[i, j] > grad_h[i, j] else grad_h[i, j]

        N_cand_mean = cv2.mean(N_cand)[0]
        mask = N_cand > N_cand_mean
        mask = mask / 255
        mask = np.float32(mask)
        N = cv2.multiply(N_cand, mask)

        sum_noise = cv2.sumElems(N)[0]
        sum_nosie_cnt = cv2.countNonZero(N)
        noise_mean = sum_noise / (sum_nosie_cnt + 0.0001)
        noise_ratio = sum_nosie_cnt / (len(N) * len(N[0]))
        return noise_mean, noise_ratio

    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    blur_mean, blur_ratio = compute_blur_IQA(image)
    # print(blur_mean,blur_ratio)
    noise_mean, noise_ratio = compute_noise_IQA(image)
    # print(noise_mean,noise_ratio)

    return 1 - (blur_mean + 0.95 * blur_ratio + 0.3 * noise_mean + 0.75 * noise_ratio)

def JPEQ(image_path):
    """

    :param image_path:
    :return:
    """

    def h_feature(x,M,N):

        d_h = np.zeros((M,N-1),np.float32)
        for j in range(N-1):
            d_h[:,j] = x[:,j+1] - x[:,j]

        DEL = math.floor(N/8) - 1

        d1_h = np.zeros((M,DEL),np.float32)

        for j in range(DEL):
            d1_h[:,j] = abs(d_h[:,8*j+7])

        B_h = cv2.mean(d1_h)[0]
        A_h = (8*cv2.mean(abs(d_h))[0]-B_h)/7

        sig_h = np.sign(d_h)

        left_sig = np.zeros((M,N-2),np.float32)
        right_sig = np.zeros((M,N-2),np.float32)

        for j in range(N-2):
            left_sig[:,j] = sig_h[:,j]
            right_sig[:,j] = sig_h[:,j+1]

        multi_sig = np.zeros((M,N-2),np.float32)

        for i in range(M):
            for j in range(N-2):
                template = left_sig[i,j] * right_sig[i,j]
                if template<0:
                    multi_sig[i,j] = 1
                else:
                    multi_sig[i,j] = 0

        Z_h = cv2.mean(multi_sig)[0]
        return A_h,B_h,Z_h

    def v_feature(x,M,N):

        d_v = np.zeros((M-1,N),np.float32)
        for i in range(M-1):
            d_v[i,:] = x[i+1,:]-x[i,:]

        DEV = math.floor(M/8)-1

        d1_v = np.zeros((DEV,N),np.float32)
        for i in range(DEV):
            d1_v[i,:] = abs(d_v[8*i+7,:])

        B_v = cv2.mean(d1_v)[0]
        A_v = (8*cv2.mean(abs(d_v))[0]-B_v)/7

        sig_v = np.sign(d_v)

        up_sig = np.zeros((M-2,N),np.float32)
        down_sig = np.zeros((M-2,N),np.float32)

        for i in range(M-2):
            up_sig[i,:] = sig_v[i,:]
            down_sig[i,:] = sig_v[i+1,:]

        multi_sig = np.zeros((M-2,N),np.float32)
        for i in range(M-2):
            for j in range(N):
                template = up_sig[i,j] * down_sig[i,j]
                if template<0:
                    multi_sig[i,j] = 1
                else:
                    multi_sig[i,j] = 0

        Z_v = cv2.mean(multi_sig)[0]

        return A_v,B_v,Z_v

    def combine_fearure(A_h,B_h,Z_h,A_v,B_v,Z_v):
        A = (A_h+A_v)/2
        B = (B_h+B_v) / 2
        Z = (Z_h+Z_v)/2

        return A,B,Z

    def compute_score(A,B,Z):
        alpha = -245.8909
        beta = 261.9373
        gamma1 = -239.8886
        gamma2 = 160.1664
        gamma3 = 64.2859

        score = alpha+beta*(pow(B,gamma1/10000)*pow(A,gamma2/10000)*pow(Z,gamma3/10000))

        return score

    if isinstance(image_path, str):
        image = preprocess(image_path)
    else:
        image = image_path
    image = np.float32(image)
    M,N = len(image),len(image[0])

    A_h, B_h, Z_h = h_feature(image,M,N)
    A_v, B_v, Z_v = v_feature(image, M, N)
    A,B,Z = combine_fearure(A_h,B_h,Z_h,A_v,B_v,Z_v)

    score = compute_score(A,B,Z)

    return score

def JPEQ2(image_path):
    """
    :param image_path:
    :return:
    """
    S = JPEQ(image_path)
    return 4.0 / (1+math.exp(-1.0217*(S-3))) + 1.0

# print(JPEQ2("../main_v1/instance.png"))
# exit()